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Abstract 

We present Car-Parrinello molecular dynamics simulations of the initial reaction steps lead- 
ing to an inert oxide layer on aluminium. The mechanism of the reaction of the aluminium 
surface with single oxygen molecules is analysed. After adsorption at the surface the oxygen 
molecules dissociate at a femtosecond timescale and the atoms are chemisorbed at the surface 
at a distance of several angstrom which falsifies earlier STM results. When the aluminium sur- 
face is exposed to higher oxygen pressure, a surface layer essentially consisting of threefold 
coordinated oxygen atoms starts to form. 
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Aluminium is among the most broadly used materials and its world production is steadily in- 
creasing.^ Its low specific weight and its resistance to corrosion render it optimally suitable for 
a large variety of applications, for example in transportation, construction and packaging. How- 
ever, as a base metal it gains its corrosion resistance from a thin layer consisting of a mixture of 
aluminium oxide and aluminium hydroxide which forms immediately when a clean aluminium 
surface is exposed to air. The composition, thickness and properties of this chemically passivating 
layer depend on the particular conditions of formation. Spontaneous self -passivation leads to lay- 
ers of a few nanometers thickness. Significantly thicker and harder layers essentially consisting of 
well-ordered a-Al203 can be obtained electrochemically. Like for bulk aluminium oxide, the sur- 
face of passivated aluminium is covered with OH groups with different coordination to the bulk.^^ 
The reaction of oxygen with aluminium was investigated with STM experiments in 1992 empha- 
sizing the role of 'hot adatoms'.f^ A subsequent study ^ demonstrated the relatively low transient 
mobility of the adsorbed oxygen atoms (0.5 nm on average). A recent HRTEM study investigates 
the growth of an aluminium oxide layer in contact with the melt.'^ 

In the present study we want to simulate, at first-principles level, the initial steps of the for- 
mation of such layers. Previous theoretical investigations using density functional theory (DFT) 
focussed on alumina surfaces and their reactivity .1^^ Early work confirmed the experimental find- 
ing that the (0001) surface is most stable .'^First-principles molecular dynamics simulations* ^ ^^^^ ^ 
showed the facile reaction of the oxide surface with water molecules leading to OH coverage. 
In a study of the reaction of oxygen with aluminium^^^^ the authors induced an artificial barrier 
to chemisorption invoking non-adiabatic effects in order to explain the low sticking coefficient 
derived from earlier experimental results. However, in view of the more recent experiments 
which indicate a much lower transient mobility and faster chemisorption,^ an artificial extension 
of DFxISini jxiay not be necessary to model the system. To analyse the mechanism we investigate 
the reaction of an aluminium surface with one oxygen molecule only simulating the gas-phase 
situation as well as the reaction with liquid oxygen. 

A series of simulations with one attacking oxygen only was performed at a temperature of 300 



2 



K using different initial orientations of the oxygen molecule relative to the surface. The incident 
oxygen molecule is moving towards the surface with a velocity of 400 m/s which leads to a reaction 
within a few hundred femtoseconds in all ten simulation runs. In each case, the oxygen molecule 
binds to the surface and dissociates. [T] illustrates the motion of the two oxygen atoms for one of 
the simulation runs. The molecule hits the surface 120 fs after the end of the equilibration. At this 
point it is accelerated to a velocity of 1300 m/s. The distance plot shows that it starts to dissociate 
but 100 fs later. Once the first oxygen atom (red graph in[TJ lower plot) contacts the surface, the Al- 
O-O angle decreases and the second atom (black graph) is bound to a threefold coordinated surface 
lattice site. During the bond dissociation the first atom is pushed over the bound aluminium atom 
to its new position while the aluminium atom itself is pulled out of the surface to another position. 
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Figure 1 : Reaction of an oxygen molecule with an aluminium surface as observed in one of the 
simulation runs. The upper plot shows the velocities of the center of mass of the two oxygen atoms 
and their distance. In the lower plot the motion of the oxygen atoms on the aluminium surface 
are depicted (black and red). One of the surface atoms (blue) is moved to another lattice site as a 
consequence of the impact. 
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[2] shows the increase in temperature of the total system. Both the adsorption and the consecu- 
tive dissociative reaction with the surface lead to a significant increase in temperature. From the 
comparison with [T] it is obvious that initially the oxygen atoms gain kinetic energy, while in the 
further course of the reaction the increase of the kinetic energy is taken up by the surface. The 
Kohn-Sham energy is lowered by roughly 0.14 a.u. (370 kJ/mol). Since no thermostats were used 
in the simulation, the electronic system heats up quite a bit and gains kinetic energy. Some exci- 
tation of electronic degrees of freedom is to be expected in this highly reactive system, however, 
since with Car-Parrinello molecular dynamics the motion of the orbitals is described classically, 
this increase in kinetic energy can be interpreted only qualitatively. The single reaction steps can 
be followed from the graphs of the charge and the spin charge Q. During the first reaction step, 
which is the adsorption to the surface, the charge of the oxygen molecule changes only partially. 
The spin charge is transferred to the surface within 160 fs while the oxygen atoms are fully ionized 
after about 400 fs. 

[3] shows some snapshots from one of the simulations. Before the start of the reaction the oxygen 
molecule is in its triplet ground state which is reflected by a high spin density. The spin density 
is transferred to the surface while the molecule binds to the surface and forms a three-membered 
ring with an aluminium atom. (The total spin of the system stays 1). After the spin charge of the 
oxygen atoms has decayed to zero, they dissociate (fourth snapshot in[3]). Some 30 fs later there is 
again a certain accumulation of spin density at the oxygen atoms which is obvious also from the 
peak in the graph of the spin charge ([2]) about 270 fs after the end of the equilibration. After this 
oscillation the oxygen atoms relax in surface lattice sites whereby an aluminium atom is strongly 
distorbed and dislocated from its lattice site to another one. 

|4] shows the variance of the O-O distances during the reaction. The motion of the oxygen 
atoms ends up at distances of about 1 a^, \/3 as, 2 a^, and a/7 as with as = 2.86 A being the 
nearest neighbor spacing. The connection to these channels should be considered with caution as 
is illustrated by the exemplary simulation run discussed above ([T]): Due to the strong disturbance of 
the surface layer, a final distance of 5 A is reached, while from the distance of the relevant lattice 
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Figure 2: Change of energy, temperature, charge and spin charge during the strongly exothermic 
reaction of a single oxygen molecule with the surface. Top graph: Kohn-Sham energy (black), clas- 
sical energy (red), total energy of the Car-Parrinello Lagrangian (green), and temperature (blue). 
The temperature is proportional to the kinetic energy of the ions and hence to the difference be- 
tween the black and the red curve, while the difference between the red and the green curve cor- 
responds to the fictitious kinetic energy of the orbitals in Car-Parrinello theory. The two lower 
graphs show the charges and spin charges of the two oxygen atoms. While the spin charge de- 
creases during adsorption, the charge of the oxygen atoms is increased during the full reaction 
including dissociation and relaxation. 
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Figure 3: Snapshots from a simulation run showing a single oxygen molecule reacting with the 
surface. The oxygen molecule (red) approaches the surface and dissociates. The spin densities are 
shown in orange. 



The average distance obtained in the ten simulation runs is 0.4 nm. This result is in nice 



adsorption of 0.5 nm. The deviation of the numerical value is within the error which stems from 
the limitations in statistics, simulation time and simulation cell size. The computed value certainly 
does not agree, however, with the value of 8 nm reported in the work by Brune et al. Obvi- 
ously the interpretation of the STM results was superficial and erroneous. The distribution of the 
adatoms was found to be random at low oxygen pressure while at higher oxygen pressure the for- 
mation of islands was observed. El Also increasing the temperature supports this island formation.^^ 
From these experimental observations and from our results, two effects may contribute to island 
formation: the initial motion of the oxygen atoms till they are chemisorbed and the strong increase 
of kinetic energy in the upper aluminium layer leading to partial melting of the metal. A motion of 
the strongly bound oxygen atoms on a 'cold' aluminium lattice seems much less likely. 

The beginning of a surface layer formation was studied in a simulation with liquid oxygen. [5] 
shows some snapshots from a molecular dynamics simulation at 300 K. After the equilibration. 



sites a distance of 7.6 A would be computed. 



agreement with the publication by Schmid et who report a mean interatomic distance after 
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Figure 4: Change of the distance of the two oxygen atoms in the ten simulations and average value. 
The average distance after dissociation and relaxation (0.4 nm) agrees well with the experimental 
value (0.5 nm). 

the distance between the nearest oxygen atoms and the surface is about 2.3 to 2.4 A. Nine out 
of 24 oxygen molecules in the simulation cell react during the first 10000 steps (484 fs), four 
of them at the slab surface shown in the figure. A thermostat was used which strongly reduces 
the kinetic energy set free in this extremely exothermic surface reaction. Nevertheless, the strong 
disorder of the top-most aluminium layer leading to partial melting at a longer timescale is obvious. 
Upon adsorption and dissociation, the oxygen atoms start to form a rigid surface layer of triply 
coordinated adatoms. 
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Figure 5: Snapshots from the molecular dynamics simulation of the reaction of the aluminium 
surface with molecular oxygen. 
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In conclusion, Car-Parrinello molecular dynamics simulations of the reaction of molecular oxy- 
gen with aluminium show that the oxygen molecules are chemisorbed immediately upon contact 
with the surface. Dissociation leads to adatoms which are separated by 0.4 nm on average in very 
good agreement with experiment. ^ Earlier experimental results are proven to be wrongl^such as the 
attempts to explain these misinterpreted experiments by artificial extensions of density functional 
theory .l^^milini Density functional theory within the BLYP approximation turns out to be very 
well suited to describe this surface reaction. The approach may serve to investigate and potentially 
falsify many more experiments in this field. 

Methods 

For our molecular dynamics simulations we used the implementation of the Car-Parrinello molecu- 
lar dynamics (CPMD) scheme!^ in the CPMD code This scheme uses density functional theory 
(j)PX)l21122| fQj. ^Y^Q description of the electronic structure. The unrestricted formulation of the 
Becke-Lee- Yang-Parr (BLYP)^^ generalized gradient approximation (GGA) functional was 
used. In the plane wave code CPMD, pseudopotentials are used for the description of the core 
electrons. The separable dual space Gaussian pseudopotentials by Goedecker, Teter and Hutter 
(GTH)™ were used with a plane-wave cutoff of 90 Rydberg. The fictitious electronic mass 
was set to the default value of 400 atomic units (a.u.). A small time step of 2 a.u. (0.048 fs) was 
chosen. The temperature of the nuclei was set to 300 K. The model system consists of four layers 
of 16 aluminium atoms stacked in an ABC A order in an orthorhombic simulation cell describing 
a (1 1 l)-surface of fee aluminium. A lattice constant of ao = 4.04959 A was used for calculating 
the cell parameters, corresponding to a fee nearest neighbour distance of as = a^^j^/l — 2.86349. 
Between the layers a spacing of approximately 12.0 A is introduced. The resulting cell parameters 
are 21.35213 A (4^^0 + 12.0), 11.45396 A (4^ao) and 9.91942 A (2y|ao). For the simulation 
of liquid oxygen, the spacing between the layers is filled with up to 24 oxygen molecules corre- 
sponding to a density of approximately 935 which is a bit lower than that of liquid oxygen 
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(1120 ^^). As liquid oxygen is highly reactive and would react immediately, it was replaced 
by unreactive nitrogen molecules during the equilibration. The temperature of the reacting sys- 
tem would rise rapidly in the liquid oxygen simulation, hence Nose-Hoover thermostats'^^^^^S'^gj-g 
used to control the temperature of the nuclei as well as the fictitious kinetic energy of the electrons. 
As thermostat parameters we use a frequency of 3000 cm~^ for coupling the nuclei to the bath and 
a frequency of 10000 cm~^ for the electrons. The fictitious kinetic energy of the electrons was 
chosen to be 0.07 a.u. The charges and spin charges were calculated by integrating the densities 
and spin densities, respectively, using Bader analysis to determine the integration range.!^ 
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